from urllib.parse import quote
from urllib import request
import json
import numpy as np
import geohash
import math
pi = np.pi
ee = 0.00669342162296594323  # 偏心率平方
a = 6378245.0  # 长半轴


def out_of_china(lng, lat):
    return not (lng > 73.66 and lng < 135.05 and lat > 3.86 and lat < 53.55)


def _transformlat(lng, lat):
    ret = -100.0 + 2.0 * lng + 3.0 * lat + 0.2 * lat * lat + \
          0.1 * lng * lat + 0.2 * math.sqrt(math.fabs(lng))
    ret += (20.0 * math.sin(6.0 * lng * pi) + 20.0 *
            math.sin(2.0 * lng * pi)) * 2.0 / 3.0
    ret += (20.0 * math.sin(lat * pi) + 40.0 *
            math.sin(lat / 3.0 * pi)) * 2.0 / 3.0
    ret += (160.0 * math.sin(lat / 12.0 * pi) + 320 *
            math.sin(lat * pi / 30.0)) * 2.0 / 3.0
    return ret


def _transformlng(lng, lat):
    ret = 300.0 + lng + 2.0 * lat + 0.1 * lng * lng + \
          0.1 * lng * lat + 0.1 * math.sqrt(math.fabs(lng))
    ret += (20.0 * math.sin(6.0 * lng * pi) + 20.0 *
            math.sin(2.0 * lng * pi)) * 2.0 / 3.0
    ret += (20.0 * math.sin(lng * pi) + 40.0 *
            math.sin(lng / 3.0 * pi)) * 2.0 / 3.0
    ret += (150.0 * math.sin(lng / 12.0 * pi) + 300.0 *
            math.sin(lng / 30.0 * pi)) * 2.0 / 3.0
    return ret


def gcj02_to_wgs84(lng, lat):
    if out_of_china(lng, lat):
        return [lng, lat]
    dlat = _transformlat(lng - 105.0, lat - 35.0)
    dlng = _transformlng(lng - 105.0, lat - 35.0)
    radlat = lat / 180.0 * pi
    magic = math.sin(radlat)
    magic = 1 - ee * magic * magic
    sqrtmagic = math.sqrt(magic)
    dlat = (dlat * 180.0) / ((a * (1 - ee)) / (magic * sqrtmagic) * pi)
    dlng = (dlng * 180.0) / (a / sqrtmagic * math.cos(radlat) * pi)
    mglat = lat + dlat
    mglng = lng + dlng
    return [lng * 2 - mglng, lat * 2 - mglat]


def isRayIntersectsSegment(poi, s_poi, e_poi):  # [x,y] [lng,lat]
    # 输入：判断点，边起点，边终点，都是[lng,lat]格式数组
    # print(poi, s_poi, e_poi)
    if s_poi[1] == e_poi[1]:  # 排除与射线平行、重合，线段首尾端点重合的情况
        return False
    if s_poi[1] > poi[1] and e_poi[1] > poi[1]:  # 线段在射线上边
        return False
    if s_poi[1] < poi[1] and e_poi[1] < poi[1]:  # 线段在射线下边
        return False
    if s_poi[0] < poi[0] and e_poi[0] < poi[0]:  # 线段在射线左边
        return False
    if s_poi[1] == poi[1] and e_poi[1] > poi[1]:  # 交点为spoint
        return False
    if e_poi[1] == poi[1] and s_poi[1] > poi[1]:  # 交点为epoint
        return False
    xseg = e_poi[0] - (e_poi[0] - s_poi[0]) * (e_poi[1] - poi[1]) / (e_poi[1] - s_poi[1])  # 求交
    if xseg < poi[0]:  # 交点在射线起点的左侧
        return False
    return True  # 排除上述情况之后


def isPoiWithinPoly(poi, epoly):
    sinsc = 0  # 交点个数
    for i in range(len(epoly)):  # [0,len-1]
        s_poi = epoly[i]
        e_poi = epoly[(i+1) % len(epoly)]
        if isRayIntersectsSegment(poi, s_poi, e_poi):
            sinsc += 1  # 有交点就加1
    return True if sinsc % 2 == 1 else False


def find_poi_id(cityname, keywords):
    amap_web_key = 'f7e5f805c9d8038065a11cf544144d94'
    poi_search_url = "http://restapi.amap.com/v3/place/text"
    i = 1
    poilist = []
    id = ''
    try:
        while True:
            page = i
            req_url = poi_search_url + "?key=" + amap_web_key + '&extensions=all&keywords=' + quote(
                keywords) + '&city=' + quote(cityname) + '&citylimit=true' + '&offset=25' + '&page=' + str(
                page) + '&output=json'
            with request.urlopen(req_url) as f:
                data = f.read()
                data = data.decode('utf-8')
            result = data
            result = json.loads(result)  # 将字符串转换为json
            if result['count'] == '0':
                break
            poilist.extend(result['pois'])
            i = i + 1
            if i > 0:
                break
        pois = poilist
        id = pois[0]['id']
    except:
        print('No id:', cityname, keywords)
        return id
    return id


def find_poi_boundary(id, cityname, keywords):
    if len(id) == 0:
        return ''
    req_url = f'https://gaode.com/service/poiInfo?query_type=IDQ&pagesize=20&pagenum=1' \
        f'&qii=true&cluster_state=5&need_utd=true' \
        f'&utd_sceneid=1000&div=PC1000&addr_poi_merge=true&is_classify=true&zoom=11&id={id}'
    f = request.urlopen(req_url)
    data = f.read()
    data = data.decode('utf-8')
    try:
        datajson = json.loads(data)  # 将字符串转换为json
        datajson = datajson['data']
        datajson = datajson['poi_list']
        datalatlng_ls = datajson[-1]['domain_list'][-1]['value']
    except:
        print('no poi boundary!!', cityname, keywords)
        return ''
    return datalatlng_ls


def gen_poi_geohash(datalatlng_ls):
    if len(datalatlng_ls) == 0:
        return ''
    else:
        arr = np.array([[float(a.split(',')[0]), float(a.split(',')[1])] for a in datalatlng_ls.split('_')])
        xmin, xmax = min(arr[:, 0]), max(arr[:, 0])
        ymin, ymax = min(arr[:, 1]), max(arr[:, 1])
        step = min(ymax-ymin, xmax-xmin)/10
        xv, yv = np.arange(xmin, xmax, step), np.arange(ymin, ymax, step)
        xv1, yv1 = np.meshgrid(xv, yv)
        xv2, yv2 = xv1.reshape(-1), yv1.reshape(-1)

        sel = []
        xy_arr = np.array([xv2, yv2]).T
        for xy in xy_arr:
            sel += [isPoiWithinPoly(xy, arr)]

        res_arr = xy_arr[sel]
        geo_ls = []
        for (lng, lat) in res_arr:
            lng84, lat84 = gcj02_to_wgs84(lng, lat)
            geo_ls += [geohash.encode(lat84, lng84)]
        geo_ls_str = ';'.join(geo_ls)
        return geo_ls_str


def proc_poi(addr, nd):
    cityname, keywords = addr
    id = find_poi_id(cityname, keywords)
    datalatlng_ls = find_poi_boundary(id, cityname, keywords)
    geo_ls = gen_poi_geohash(datalatlng_ls)
    if len(geo_ls) > 0:
        geo_ls1 = geo_ls.split(';')
        geo_ls2 = ';'.join(list(set([a[:nd] for a in geo_ls1])))
    else:
        return ''
    return geo_ls2


if __name__ == '__main__':
    # cityname = "湖北省"
    # keywords = "湖北省"
    # nd = 8
    # geo_ls = proc_poi((cityname, keywords), nd)

    import pandas as pd
    import sys

    if len(sys.argv) < 2:
        file_in = 'data/yiyuan.xlsx'
        file_out = 'res/yiyuan_res.xlsx'
        nd = 8
    else:
        file_in = sys.argv[1]
        file_out = sys.argv[2]
        nd = int(sys.argv[3])

    A = pd.read_excel(file_in, header=0)
    A['geo_ls'] = A.apply(proc_poi, axis=1, **{'nd': nd})
    A.to_excel(file_out, header=True, index=False)
    print('输入数据量: ', A.shape[0], '缺失poi数据量: ', sum(A['geo_ls'] == ''))


